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ABSTRACT 

We examine the jets and the disc of SS 433 at super-Eddington luminosities with 
M ~ 600Mc by time-dependent two-dimensional radiation hydrodynamical calcula- 
tions, assuming a-model for the viscosity. One-dimensional supercritical accretion disc 
models with mass loss or advection are used as the initial configurations of the disc. 
As a result, from the initial advective disc models with a =0.001 and 0.1, we obtain 
the total luminosities ^ 2.5 x 10"*° and 2.0 x 10*° erg s^^. The total mass-outflow rates 
are ~ 4 x 10^^ and 1O^'*M0 yr^^ and the rates of the relativistic axial outflows in a 
small half opening angle of ~ 1° are about 1O~^M0 yr~^: the values generally con- 
sistent with the corresponding observed rates of the wind and the jets, respectively. 
From the initial models with mass loss but without advection, we obtain the total 
mass-outflow and axial outflow rates smaller than or comparable to the observed rates 
of the wind and the jets respectively, depending on a. In the advective disc model 
with a = 0.1, the initially radiation-pressure dominant, optically thick disc evolves to 
the gas-pressure dominated, optically thin state in the inner region of the disc, and 
the inner disc is unstable. Consequently, we find remarkable modulations of the disc 
luminosity and the accretion rate through the inner edge. These modulations manifest 
themselves as the recurrent hot blobs with high temperatures and low densities at the 
disc plane, which develop outward and upward and produce the QPOs-like variability 
of the total luminosity with an amplitude of a factor of ^ 2 and quasi-periods of ^ 10 
- 25 s. This may explain the massive jet ejection and the QPOs phenomena observed 
in SS 433. 

Key words: accretion, accretion discs - black hole physics - hydrodynamics - radi- 
ation mechanism: thermal- X-rays:individual: SS 433. 



1 INTRODUCTION 

Disc accretion is an essential process for such phenomena 
as energetic X-ray sources, acti ve galactic nucle i , and pro- 
tostars. Since the ea r ly wor ks bv lPringle fc Ree3 (|l972l ') and 
IShakura fc SunvaevI (l973l ), a great number of papers have 
been devoted to studies of the disc accretion onto gravi- 
tating objects. When the accretion rate is not too high, 
the accretion disc luminosity is directly in proportion to 
the accretion rate and can be successfully described by the 
Shakura-Sunyaev (S-S) model. However, for the supercriti- 
cal accretion rate, matter flows out of the disc, and the rate 
of accretion onto the central black hole is reduced, regu- 
lating the luminosity to the Eddington limit. The galactic 
microquasar SS 433 is a promising super-critically accret- 
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ing sterllar-mass black hole candidate and has stimulated 
numerous studies, because it displays remarkable observa- 
tional features, such as its extremely high energy, two op- 
positely directed relativistic jets, and the precessing mo- 
tion of the jets. Although a number of observational and 
theoretical studies on SS 433 have been published (for de- 
tailed reviews see Margon 1984; Fabrika 2004), there are 
still many problems of the jets and the disc to be solved. 
The super-Eddington accretion discs are generally expected 
to possess vortex funnels and radia tion-pressure driven 
jets fr om geometr i cally t hick discs ('Shakura fc Sunyaej 
19731: iLvnden-Belll 1 19781 : iFukud [l982 : Calvani fc Nobiii| 



19831 ) . Such supercritical disc models have been numerically 



examined by two-dimensional ra diation hy drodynamical 
calcu l ations (lEggum. Coroniti fc~k atz 19851. Il988l: lOkudal 



2002'; Okuda ot al."2005' jOhsuga et al.l 120051 : lOhsugal 12007 
Ohsuga fc Mineshigai2007h . especially focussing on the puz- 
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Figure 1. Solid lines: accretion rate m(r) normalized to the input 
accretion rate and relative thickness h/r oi the supercritical disc 
with mass loss for the case of the input accretion rate mo ~ 
600 and the viscosity parameter a = 0.1. Dashed lines show the 
solution for S-S model without mass loss. 
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Figure 2. Same as Fig.l but for the profiles of temperature T 
(K) and density p (g cm~^). The temperature profile is the same 
for the both models, as far as the radiation-pressure dominant 
disc is concerned. 



zling X-ray sourc e SS 433 | Eggum, Coroniti fc Katdll985l . 
Ll988; Okuda 2002: IOkudaet al. 2005). However, these works 
leave something to be desired as far as SS 433 is concerned, 
because the accretion rates adopted in these studies are very 
small, compared with those estimated for SS 433. In this pa- 
per, we examine the properties of the jets and the disc of 
SS 433 with the plausible supercritical accretion rate. We 
perform time-dependent two-dimensional radiation hydro- 
dynamical calculations, using ID models of the su percritical 
accre tion disc models with mass loss or advection l|Lipunoval 
1 19991 ) as initial disc configurations. 



Figure 3. Disc thickness h/r, temperature T (K), and density p 
(g cm^'^) of the supercritical disc with advection for the case of 
mo ^ 600 and a = 0.1. 

2 SUPERCRITICAL ACCRETION DISC 
MODELS WITH MASS LOSS AND 
ADVECTION 



Our previous studies of SS 433 (|Okudall200i lOkuda et al.l 
l2005h are based on the initial discs by S-S model with 
rho ^ 20. Here, mo is the input accretion rate normalized 
to the Eddington critical accretion rate Mc(= L^/r)(? = 
A-kGM /r]CK), where Le is the Eddington luminosity, G the 
gravitational constant, M the black hole mass, k the Tho- 
moson scattering opacity, c the velocity of light, and rj the 
conversion efficiency of the energy of the accreting matter 
into radiation. Adopting 77 — 1/12, we have Mc = 2.0 x 10^^ 



g s 



and Le = 1.5 x 10' 



,39 



erg s"' for M = IOMq. In 
600, which corresponds to 
"^Mp) yr~ "^ and is a more plausible value 
Fabrika 2004'). 



this paper we consider mo 
Mo ^ 1.8 X 10 
for SS 433 



The S-S disc with such a high 
accretion rate is geometrically too thick and becomes in- 
valid against assumptions used. When the accretion rate is 
very high, the disc luminosity exceeds the Eddington lu- 
minosity and matter flows out of the disc surface, causing 
the accretion rate onto the central object to decrease. The 
outflow takes place inside the radius at which the disc thick- 
ness becomes comparable to the disc radius. This radius is 
called a sp herization radius Rsp (jShakura &: Sunvaevlll973l ). 
iLipunoval ^1999) proposed supercritical disc models, taking 
account of mass loss and advection through the accretion 
discs, which were applied to SS 433 and the ultr aluminous 
X-ray sources (ULXs) by iPoutanen et al] l|2007h . In their 
models with mass loss, it is assumed that a fraction ew of 
the radiation energy flux is spent on the production of the 
outflow within the spherization radius. For ew = 1 , if advec- 
tion is neglected, we obtain analytical solutions l|Lipunoval 
I1999I ) for the disc variables at the disc plane, and the ra- 
tio Rsp/3Rg is estimated to be ~ |mo, where Eg is the 
Schwarzschild radius. The accretion rate in the disc within 
the spherization radius is given by 

x-5/2 



M(r) = Mo 



1 + 



2 ( flap 

3 '> 3R, 



-5/2 



for r < Rs 



(1) 



This law has the asymptotic form 
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Table 1. Results and comparison with observations 
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Figure 4. Time evolution of total luminosity L, disc luminosity L^, total mass-outflow rate Mout ejected from the entire outer boundary 
surface, mass-outflow rate Mio of the axial outflow with a half opening angle 1.2° along the rotational axis, and mass-inflow rate Mcdgc 
into the black hole through the inner boundary for model ML-2, where time is shown in units of Rg/c. 



Mir) « Mo for 3i?g < r i?ep. (2) 

The temperature T and the density p in the central plane 
of the disc are approximately given in the region of 3i?g ^ 
r < Rap by 

-<^'— '(T^)""(^)-"fe)""''<'"> 

where a is the viscosity parameter. 

Figs 1 and 2 show the accretion rate m(r), the relative 
disc thickness h/r, the central plane temperature T, and the 
density p for the ID model with mass loss and the viscosity 
parameter a — 0.1, along with the solution for S-S model 
without mass loss. Outside the spherization radius the solu- 
tions are identical. 

The disc at the supercritical accretion rates be- 
comes geometrically thick. The emission from the sur- 
face of a thick disc is not an efficient cooling mecha- 
nism, and the advective transport of the viscously gener- 
ated heat becomes important in the energy balance equa- 
tion (|Paczvriski fc Bisnovatvi-KoganI ,1981). An advection- 
dominated disc model at moderately supe r-Eddington accre- 
tion r ates mo 5; 70 was first proposed bv lAbramowicz et al.l 
j 19881 ). A supercriti cal disc wi t h even higher accretion rates 
was approached by iLipunoval l|l999l ). where mass loss and 



advection were included, and the discs were geometrically 
thick [h/r ~ 0.7 at maximum) but not so thick as it followed 
from S-S model (h/r 100). Fig. 3 shows a solution for the 
advective supercritical disc with mass loss for rho ~ 600 
and Q = 0.1. This solution and the analytic solution for the 
disc with mass loss, described above, are used in the current 
work as initial configurations of the disc in the 2D numerical 
model. 



A set of relevant equations for the numerical calculation con- 
sists of six partial differential equations for the density, the 
momentum, and the thermal and radiation energy. These 
equations include the full viscous stress tensor, the heat- 
ing and cooling of the gas, and the radiation transport . 
The pseudo-Newtonian potential (|Paczyriski fc Wiitalll980h 
is adopted in the momentum equation. The radiation trans- 
port is t reated in the gray, flux-limited diffusion approxi- 
mation l|Levermore fc Pomranin^ Il98ll ). We use spherical 
polar coordinates {r,^,ip), where r is the radial distance, C 
is the angular distance measured from the equatorial plane 
of the disc, and (p is the azimuthal angle. The above set 
of differential equations is numerically integrated in time 
using a finite-difference method, wh ich is an improved ver- 
sion of that described in iKlevI (Il989h . The method is based 
on an explicit-impli cit finite difference scheme, whose de- 
tails are described bv lOkuda. Fuiita fc Sakashital |l993) and 



\ -1/2 

^1 gcm~^(4) 
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lOkuda et al] l|2005h . The computational domain is divided 
into Nr X A^,; grid cells, where Nr grid points (=100) in the 
radial direction are spaced logarithmically as Ar/r = 0.1 
and A^,; grid points (=150) in the angular direction are 
equally spaced, but more refined near the equatorial plane, 
typically with AC, = 7r/150 for 7r/2 ^ C ^ '"'/S and 
A(" = 7r/300 for 7r/6 5^ ^ 0, in order to resolve the struc- 
ture of the accretion disc. Although the radial mesh-sizes do 
not have a fine resolution to examine detailed disc structure, 
we consider the mesh-sizes to be sufficient for examination 
of the global behavior of the disc, the jets, and the wind 
mass-outflow rates. 



3.1 Model Parameters 

We consider a Schwarzschild black hole with mass M — 
IOM0 and take the inner-boundary radius i?in = 2J?g and 
the outer boundary radius of the spherical computational 
domain Rout = 4.5 x 10''-Rg(~ 1.4 x 10^^ cm). For the kine- 
matic viscosity, we adopt the usual a - model. In Table 1, 
we give the parameters of the discs with mass loss (ML) 
and advection (AD), where Mo is the input accretion rate, 
and a is the viscosity parameter. Starting with ID solutions 
described in section 2, we perform time integration of the 
equations until a quasi-steady solution is obtained. 



3.2 Initial Conditions 

The initial conditions consist of a dense, optically thick ac- 
cretion disc and a rarefied optically thin atmosphere around 
the disc. Physical variables p and T at the equatorial plane 
for r ^ 3-Rg are given by the ID solutions for the super- 
critical discs with mass loss or advection. We construct the 
vertical structure of the disc in the approximation of the 
hydrostatic and radiative equilibrium by integrating the rel- 
evant equations with given boundary values at the equato- 
rial plane. As for the initial atmosphere around the disc, 
assuming that the gas is in the optically thin limit and in 
the radially hydrostatic equilibrium, we have 



Ac dEr 
pKc dr 



dP 



GM 

2~' 



(5) 

(6) 



or r 

where Fr is the radial component of the radiative flux F, 
Er the radiation energy density, A the flux-limiter, Kc the 
electron-scattering opacity, and P the total pressure. Fur- 
thermore, if it is assumed that the flux-limiter A is constant 
and that the radiation pressure is dominant throughout the 
gas, we have for the initial radiation energy density Ei and 
the density p, from the above equations. 



Er 



1 GM 



2A 



(7) 



(8) 



Actually, the flux-limiter A is taken to be ~ 10~^. How- 
ever, we note that a particular initial distribution of the gas 
around the disc does not influence the results at the suffi- 
ciently late simulation times. 



t-14t 
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r/Rg 



Figure 5. Radial profiles of density p (g cm~^), temperature T 
(K), radiation temperature Tr (K) (dashed line), @ and ratio /3 
of the gas pressure to the total pressure on the equatorial plane 
at t = 5.2 X lO^Rg/c for model ML-2. 



3.3 Boundary Conditions 

Physical variables at the inner boundary, except for the ve- 
locities, are calculated by extrapolation of the variables near 
the boundary. We impose limiting conditions that the radial 
velocities at the inner boundary are given by a fixed free-fall 
velocity and the angular velocities are zero. On the rota- 
tional axis and the equatorial plane, the meridional tangen- 
tial velocity w is zero and all scalar variables must be sym- 
metric relative to the axis and the plane. The outer bound- 
ary at r = Rout is divided into two parts. One is the disc 
boundary through which matter is entering from the outer 
disc. At this outer-disc boundary we assume a continuous 
inflow of matter with a constant accretion rate Mq. The 
other part is the outer boundary region above the accretion 
disc. We impose free-floating conditions on this outer bound- 
ary and allow for outflow of matter, whereas any inflow is 
prohibited here. We also assume the outer boundary region 
above the disc to be in the optically-thin limit, |i^| ^ cE^. 



4 NUMERICAL RESULTS 

In order to obtain a reliable configuration of the jet, the sim- 
ulation time should be longer than the jet arrival time tjot at 
the outer boundary. On the other hand, to get a steady state 
of a viscous disc in a region of radius r, we need the sim- 
ulation time longer than the characteristic dynamical time 
tdyn and the viscous time-scale t^is- These times are given 
as follows: 



tu 



tdy 



and 



Rou 

'vT 



2.2 x 10' 



jot 



(VjctV 
\0.2cJ 



'?(7)"'--»'(^)l^) 



3/2 



Rrr 



- 3x10* 



\^103i?g 



3/2 



Rn 



(9) 



,(10) 
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where n and Vjet are the Keplerian angular velocity and the 
typical jet velocity. For the supercritical disc of SS 433, we 
have h/r ~ 1 and Vjet ~ 0.27c. Large-scale mass loss sets in 



near the spherization radius of the disc at 



For 

this radius and a = 0.1, the maximum time of tjct, tdyn, and 
tvis is ~ 3 X 10^-Rg/c. We set this time as an approximate 
computational goal, exceeding it in the cases of models ML-2 
and AD-2 in order to examine the disc instability. Whereas 
the quasi-steady values of the luminosities are almost at- 
tained at the final phases, the simulation times are still not 
sufficient for the disc and the outflow to settle into a com- 
pletely steady state, because the input accretion rate is not 
equal to the total mass-outflow rate plus the mass-flux rate 
of the falling gas into the black hole. In spite of the lim- 
ited computational time, we are able to derive the general 
properties of the disc and the outflow. 

The luminosity curve is a good measure to check if a 
steady state of the disk and the outflow is attained. The 
total luminosity L and the disc luminosity Ld are calculated 
as J FdS, where the surface integral is taken over the outer 
boundary surface and the disc surface, respectively. The disc 
surface is placed at the height where the density drops to a 
tenth part of the value at the equatorial plane. This lead to 
some uncertainty in Ld because the vertical structure of a 
geometrically thick disc is treated rather approximately. 

To compare with the observational data for SS 433, 
we calculte the total mass-outflow rate Mout and the "ax- 
ial outflow" rate Mi" . The value Mout is calculated for the 
entire outer-boundary surface and corresponds to the ob- 
served mass-outflow in the wind of SS 433. Observations 
indicate that the opening angl e of the X-ray and optical 
jets i s ~ 1.2° ( MarshaU, Canizares fc Schuljf20o3 : iFabrikal 
|2004| ). Throughout the paper, we call the relativistic out- 
flow through the outer base of the cone directed along the 
rotational a:xis with the half opening angle of 1.2° "the ax- 
ial outflow". The mass-flux rate Afio of the axial outflow 
should correspond to the mass-outflow rate in the jets of SS 
433, provided the size of the present computational domain 
is comparable to the distance of the observed jets from the 
central source. In order to check the total mass-flux rate, 
we also calculate the mass-flux rate Mcdgo of the gas falling 
through the inner boundary into the black hole. 

In Table 1 we give L, Ld, Afout, Afio, and Afcdgo at 
the final simulation time t^v for each of the models. Val- 
ues L, Ld, and A/odgo show QPOs-like features in model 
AD-2; Table 1 shows their averaged values around the fi- 
nal phase. In the columns for Aifout and M-i_o , one can see 
the observational mass-outflow rate in the wind and the jets 

' 1 } ! 1 ' 

of SS 433 (Dopita & Cherepashchuk 1981; van den Heuve 



19811: iKotani 1998; Marshall, Canizares fc Schulz. .200 



Cherepashchuk et all 1200^ : IFabrikal l20oi ) 



S 



4.1 Cases of ML-1 and ML-2 with mass loss 

In the model ML-1, the initial disc is radiation-pressure 
dominant and optically thick. The perturbed density and 
temperature waves are generated initially near the inner 
edge, and then heating and cooling waves move up and down 
between the inner edge and the radius of ~ 8i?g, raising 
the disc temperature to ~ 10^^ K. As a result, the gas- 
pressure dominat, optically thin region is formed in the in- 
nermost region. The process is repeated for a while, then 



1500 




Figure 6. Velocity vectors and density contours near the spher- 
ization radius on the meridional plane at t = 3.2 X W^Rg/c for 
model ML-2. The reference vector of light speed is shown by a 
long arrow. The axial outflow propagate with relativistic veloci- 
ties ~ 0.2c along the Z-axis. 
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Figure 7. Same as Fig. 6 but for temperature contours. The 
velocity vectors are indicated by unit vectors. 



stops after t ~ 4 x 10''iig/c. Finally, we get stable, opti- 
cally thick, and radiation-pressure dominant disc. The re- 
sultant disc temperatures at the equatorial plane are not 
different signiflcantly from the initial temperatures of the 
disc, while the densities are enhanced by a factor of 10 - 
30 in the region of r ^ 10^ -Rg. The flnal luminosities L and 
Ld ~ 2 X 10*" erg s^^ are one order of magnitude larger 
than the Eddington luminosity. The total mass-outflow rate 
Afout ~ 4.5 X 1O~'^M0 yr"'^ is less by an order than the mass- 
outflow rate of the wind observed in SS 433, and the rate of 
the axial outflow Mio ~ 1.6 x 1O~^M0 yr~^ is marginally 
in the range of the observed mass-outflow rate of the jets. 

Fig. 4 shows the time evolutions of L, Ld, Mout, Mio , 
and Afedge for model ML-2 with the viscosity parameter 
a=0.1. The total luminosity L becomes comparable to the 
disc luminosity Ld at the time t ~ Rout/c = 4.5 x 10'*7?g/c, 
which is the photon transit time from the inner edge to the 
outer boundary. The total luminosity L 10*° erg s"'^ ob- 
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Figure 8. Time evolutions of L, Lj, Mout, M^o, and M^dgc for 
model AD-1. 



tained finally is greater by a factor of 6 than the Eddington 
luminosity. 

The initial disc of ML-2 is radiation-pressure dominant 
and optically thick throughout the whole disc region. In the 
same way as in the case of ML-1, the perturbations of den- 
sity and temperature, generated initially near the inner edge, 
propagate outward and inward as heating and cooUng waves, 
and the high temperature, gas-pressure dominant, optically 
thin region is formed in the inner disc. The instability in 
the optically thin region is never damped, as indicated by 
MedgG in Fig. 4, and appears as an oscillating hot blob with 
a variable size of 10 - 100 -Rg at the equatorial plane. After 
a large-scale hot blob at f ~ 1.6 x 10''i?g/c, the instabil- 
ity persists but the absolute variation amplitude of Mcdgc 
becomes smaller. In spite of the considerable variability of 
Medge, the modulations of L and are weak, except for 
the phase of the large-scale blob, and become negligible at 
the later phase. At the final phase, the first outward heat- 
ing wave reaches the spherization radius and merges into 
the outer Shakura-Sunyaev disc, and the whole disc tends 
to settle gradually into a steady state. Eventually, the disc 
evolves into two zones: the gas-pressure dominant and op- 
tically thin inner disc and the radiation-pressure dominant 
and optically thick outer disc. 

Fig. 5 shows the radial profiles of density p (g cm~^), 
gas temperature T (K), radiation temperature Tr (K), and 
ratio /? of the gas pressure to the total pressure on the equa- 
torial plane at t = 5.2 x 10''7?g/c for model ML-2. The ra- 
diation temperature T-c is defined as Tr — (E-c j a)^^'^ , where 
a is the radiation density constant. When the gas is fully 
optically thick, the gas temperature T should be equal to 
the radiation temperature Tr. When the gas becomes opti- 
cally thin, Tr is much lower than T as it is found in Fig. 5. 
Initially, the radiation temperature Tr and the gas temper- 
ature T are ~ lO'' K near the inner edge of the optically 
thick disc. However, after some time, the gas temperature 
near the inner edge goes up to 10^° - 10^^ K and is variable 
by a factor of 10. Similarly, the radiation temperature is also 
modulated between 10® - 10^ K in the region of r ^ 100-Rg. 
The gas-pressure dominant, hot, optically thin inner disc is 
separated from the radiation-pressure dominant outer disc 



by the transition region at 



100i?g. The gas variables 



T and p vary sharply across the transition region, but the 
radiation temperature Tr is smooth everywhere. 

Figs 6 and 7 show the contours of density p (g cm"'' ) 
and temperature T (K) with the velocity vectors of the gas 
near the spherization radius ?A,t — 3.2 x lO^iig/c for model 
ML-2. Here we see a rarefied, hot, and optically thin high- 
velocity outfiow region and a dense, cold, and geometrically 
thick disc. The outfiow propagates with relativistic velocities 
of 0.08 - 0.2c along the Z-axis. The disc is geometrically thick 
with hjr ~ 0.7 - 1 in the region of 10^i?g < r ^ 4 x \{?R^ 
and the large-scale outfiow sets in the region. In the high- 
velocity region along the rotational axis, the temperature is 
as high as ~ 10^ K and the density is as low as ~ 10~^^ - 
lO"'' g cm-^ 



4.2 Cases of AD-1 and AD-2 with advection 

Fig. 8 shows the time evolutions of L, Ld, Mout, Mio , and 
Modgc for model AD-1. The evolutionary features of the disc 
is very similar to the case of ML-1. Finally, we get stable, 
optically thick, and radiation-pressure dominant disc with 
the maximum disc temperature ~ 10* K in the innermost 
region. The total luminosities L ~ 2.5 x lO"*" and 2.0 x lO"" 
erg sT^ in models AD-1 and AD-2 are in the same range 
as those in models ML-1 and ML-2. The total mass-outfiow 
rates are ~ 3.8x 10~^ and l.Sx 10~'*Mq yr~^ and the rates of 
the relativistic axial outfiow are 9.2x10"'' and 1.2xlO"®M0 
yr"^ in models AD-1 and AD-2 respectively. These values 
of Mout and M\° are consistent with the corresponding ob- 
served rates of the wind and the jets. 

Fig. 9 shows the time evolutions of L, Ld, Mout, Mio , 
and Afcdgo for model AD-2. In this case, we find the remark- 
able variabilities of Mcdgc, Ld, and L and, to examine their 
properties, we took a long simulation time t ~ 10® -Rg/c (~ 
100 s) . The high- velocity jets propagate vertically to the disc 
plane and expand gradually from the rotational axis with in- 
creasing time. After the time i?out/0.2c, the jets arrive at the 
outer boundary in the polar direction, and the mass outfiow 
begins and gradually approaches the steady state. In model 
AD-2, most of the accreting matter is fiown out as wind and 
only 10 percent of the input matter is swallowed into the 
black hole. 

The modulations of Mcdgc in AD-2 show two types of 
variability: (1) the small amplitude variations with a short 
time-scale, (2) the large amplitude variations with a long 
time-scale. The disc luminosity is also modulated by a fac- 
tor of a few to ten, synchronously with A/cdgc. The varia- 
tions with the small amplitude and short periods strongly 
infiuence the disc luminosity but not the total luminosity 
measured at the distant outer boundary. If the atmosphere 
between the disc surface and the outer boundary was fully 
optically thin, the total luminosity would suffer from the 
same modulations as the disc luminosity. However, in the su- 
percritical accretion flow with high input density, only the 
large modulations with the long periods contribute to the 
total luminosity, due to the atmospheric absorption around 
the disc. These modulations behave as QPOs-like variabili- 
ties. In Fig. 10, we plot the power density spectra of Medge, 
Ld, and L and recognize the QPOs-like features of these 
spectra with characteristic signals of ~ 4 x lO"'^ - 10~^ 
and 0.5 - 2 Hz for Afcdgc and Ld. Only longer-period vari- 
ations are obtained for L; therefore, we expect modulations 
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Figure 9. Time evolutions of L, Lj, Mout, Mio, and Medge for model AD-2. The disc luminosity and the mass accretion rate at the 
inner edge show two types of variability: one with the small amplitude and short periods and another with large amplitude and long 
periods. Only the longer-term variations influence the total luminosity in the form of small modulations measured at the distant outer 
boundary. 
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Figure 10. Power spectra of time variations of the total lumi- 
nosity L, the disc luminosity Lj, and the mass-inflow rate M^^g^ 
at the inner edge of the computational domain for model AD-2. 



of the observed luminosity with the quasi-periods between 
~ 10 and 25 s. 

Fig. 11 shows the radial profiles of density p (g cm"'') 
and temperature T (K) at the equatorial plane for t = 
1.03 X 10^ (dashed lines) and 1.06 x lO'^RJc (solid lines) 
for model AD-2. Similarly to the case of model ML-2, the 
disc obtained in model AD-2 consists of three regions; (1) 
the gas-pressure dominant, optically thin disc in the in- 
ner region, (2) the advection-dominated, radiation-pressure 
dominant, optically thick disc in the intermediate region, 
(3) the outer radiation-pressure dominant, optically thick 
Shakura-Sunyaev disc at r ^ 10*7?g. The regions (1) and 
(2) are sharply separated by the critical radius re- The criti- 
cal radius moves randomly up and down, and usually reaches 
~ 60iig with a heating wave, advancing sometimes to ~ 100 
- 200 Rg and rarely to the maximum radius ~ 500-Rg, then 
recedes to the inner edge with a cooling wave. These insta- 
bilities in the inner disc lead to recurrent hot blobs with high 



Figure 11. Radial profiles of density p (g cm~^) and tempera- 
ture T (K) at the equatorial plane for AD-2, where dashed and 
solid lines show the profiles at the phases of t = 1.03 X 10® and 
1.06 X 10*^ Rg/c. The gas-pressure dominant, optically thin disc 
in the inner region is separated from the advection-dominated, 
radiation-pressure dominant, optically thick disc by the critical 
radius rc, which moves randomly up and down, usually reaches 
~ 60-Rg, and sometimes advances to ~ 100 — 200 ijg. 



temperatures and low densities, like bubbles in the boilling 
water. The hot blobs typically grow to the size of ~ 60-Rg in 
the disc plane, go up with increasing size, and finally decay 
at z ~ 2000-Rg . The time evolution of the hot blobs is shown 
in Fig. 12. 

Figs 13 and 14 show the contours of gas temperature T 
(K) and radiation temperature Tr (K) with velocity vectors 
in the whole computational domain at t = 1.1 x 10^ Rg/c for 
model AD-2. The gas temperature T in the high-velocity 
jets region is as high as ~ 10* K at r ^ 10^° cm and T 3> 
Tr because the jet region is mildly optically thin, but T ~ 
Tr in the outer optically thick disc. From Fig. 14, we find 
that the contours of Tr, as well as the contours of radiation 
energy density Er, show an anisotropic distribution of radial 
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Figure 12. Time evolution of hot blobs in model AD-2, where the temperature contours of the blobs are shown at t = 1.74 X 10^, 
1.76 X 10^, 1.79 X 10^, and 1.81 x lO^Rg/c. The hot blobs typically grow to the size of ~ 60-Rg in the disc plane, go up with increasing 
size, and finally decay at z ~ 2000-Rg. 



component of the radiative flux F, where F oc — grad i?r. 
Actually, the radial components -Fr in the direction oi C, ^ 
70° exceed those in the direction of ~ 30° by a factor of 5 
- 7. 

Fig. 15 shows the contours of density p (g cm~"^) with 
magnified velocity vectors over the whole computational do- 
main, where the velocity vector of 0.2c is denoted in the 
legend. From detailed analyses of the density contours and 
the velocity vectors, we recognize three characteristic out- 
flows originating in the different regions of the disc shown in 
Fig. 11: (1) the most relativistic axial outflow with velocities 
~ 0.1 - 0.3c ejected perpendicularly to the innermost hot, 
optically thin disc, (2) the high-velocity (0.1 - 0.05c) out- 
flow within a half opening angle ~ 30° ejected from the 
advection-dominated, optically thick disk in the interme- 
diate region, (3) the slow outflow with velocities ~ 0.01c 
flowing from the disc region near the spherization radius. 
The slow outflow (3) from the outer disc, interacting with 
the high- velocity outflow (2), is blown obliquely beneath the 
high- velocity flow and is accelerated up to the velocities of 
0.002 - 0.05c for 10° < C 60° at the outer boundary. 
The present result for the supercritical disc with rho ~ 600 
shows a broader opening angle of the entire outflow from 
the disc, comparing with ^ 30° in the previous study with 
rho ~ 20 lOkuda et al.ll2005l ). This is due to the reason that 
the supercritical disc with a very high accretion rate has a 
very large spherization radius, where a massive outflow (3) 
is going out distortedly from the disc. 

5 INSTABILITY OF 

ADVECTION-DOMINATED DISC 

The initial discs in all cases considered here are radiation- 
pressure dominant and optically thick throughout the whole 



disc region. Therefore, their stability may be supposedly in- 
terpreted in terms of the slim disc model at highly super- 
Eddington luminosity; that is, these discs should be sta- 
ble against the local and global perturbations. Actually, the 
time evolution of the disc in models ML-1 and AD-1 with 
small ct — 0.001 eventually shows stable features. On the 
other hand, just at the beginning of the simulations in ML- 
2 and AD-2, the unstable behaviors of the accretion rate 
at the inner edge and the disc luminosity are revealed. This 
can be attributed to the fact that the gas-pressure dominant, 
optically thin, high temperature state is triggered in the in- 
nermost region of the disc due to initial perturbations of the 
disc variables. The optically thin disc region, as that devel- 
oping repeatedly in model AD-2, is not attained finally in 
models ML-1 and AD-1, because the disc densities in these 
models are much higher than those in ML-2 and AD-2 due 
to /9 (X in the initial disc. 

The thermal instability of advection- 
dominated one-temperatu re dis cs was examined by 
iKato. Abramowicz fc ChenI l|l996l ). It was shown that, 
in the case of a - viscosity, the optically thin advective 
disc is unstable against local perturbations if the viscosity 
parameter a is small, but that two-dimensional analysis is 
necessary to investigate stability if a is large. The stability 
of an optically thin, advection-dominated accretion disc 
with large a = 0.3 wa s examined by ID tim e-dependent 
numerical simulations l|Manmoto et al. I I1996D . and it was 
found that any disturbance added onto the accretion flow 
at large radii does not decay rapidly and continue to be 
present as fluctuations in the X-ray emission of an accretion 
disc, though the global disc structure is not modified. From 
the present simulations, we find that the inner, optically 
thin, gas-pressutre dominant disc in model AD-2 becomes 
unstable eventually; the disc luminosity fluctuates, though 
the total luminosity at the outer boundary does not vary 
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Figure 14. Same as Fig. 13 but for radiation temperature Tr (K). 
If the gas is fully optically thick, Tr is equal to the gas temperature 
T (K). 



significantly. Thus, our two-dimensional calculations of the 
advection-dominated discs confirm the above theoretical 
and ID numerical results. On the other hand, for model 
ML-2 with large q, we find that i\'fcdgc, T, and Tr are 
modulated by a factor of 3 - 10 in the innermost region of 
the optically thin disc, that is, the disc is locally unstable. 
However, no significant modulation of the disc luminosity 
is found especially at the later phases. This is attributed 
to the fact that, in model ML-2, the absolute value of the 
local accretion rate in the inner region and its variation 
amplitude are small compared with the input accretion 
rate. Therefore, the modulations do not influence the lu- 
minosities of L and Ld, because the resultant gravitational 
energy release due to the accreting gas could not contribute 
largely to the total radiation. 



SS 433 is a typical stellar-mass black-hole candidate with a 
highly super-Eddington luminosity. The present results can 
be compared with the observations of SS 433, because the 
model input accretion rate Mo ( ^ 1-8 x 1O~'*M0 yr~'^) is 



in the observed range of SS 433 mass transfer rate I Fabrikal 

The unique jet velocity Vjet (=0.26 c) of SS 433 is ex- 
plained in terms of the relativistic velocity of plasma accel- 
erated by the radiation-pressure force in the inner disc. The 
X-ray spectral lines from the gas moving with relativistic 
velocities are emitted in the very hot, optically thin region 
along the rotational axis. The temperature T ~ 10^ - 10^ K 
and the density p ~ 10~* - 7x 10"^" g cm"'' obtained in the 
present simulations at distances 6 x lO" cm ^5 r ^ 1.4 x 10^^ 
cm in the conical axial outfiow are generally consistent with 
the values T ~ 1.1 x 10** - 6 x 10** K and p ~ 3 x 10"^ - 
7 X 10~^^ g cm"'', obtained by fitting the model of X-ray 
jet emission lines produced at distances up to 2 x 10 ^^ cm 
from the jet base (|Marshall. Canizares fc Schulj|2002l ). 

In models AD-1 and AD-2, we obtain the total mass- 
outflow rates Mout (3.8 x 10"^ and 1.3 x lO'^M© yr"^) 
and the axial mass-outflow rates Afio (9.2 x 10~^ and 
1.2 X 10~^Mq yr~^) comparable with the observed wind 
(~ 10-5 - IQ-'Mo yr"^) and jet (- lO"'' - W'^Mq 
yr"'^) mass-ouflow rates, respectively. On the other hand, 
the mass-outflow rate Mout in ML-1 and the axial outflow 
rate Mio in ML-2 are by more than an order of magnitude 
lower compared with the corresponding observed values of 
the wind and jets, respectively. 

The observed half opening angle of SS 433 jets is found 
to be very small, about ~ 1° , in the distant region of ~ 10^" 



lO' ^^ cm for the X-ray jets (|Marshall. Canizares fc Schulj 
120021 ) and 10^* - lO^'^ cm for optical jets l|Fabrikall2004 ). Our 



computational domain (r ~ 1.4 x 10 cm ) is located roughly 
at the base of the observed X-ray jet, and the obtained half 
opening angle 8c for the high-velocity outflow with radial 
velocities of K ^ 0.05c is rather large, 6c ^ 30°. The rela- 
tively large half opening angles of high-velocity flows, ~ 30° 
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- 60° , have been also found i n other numerical studies of the 
supercritical accret i on discs jEggum. Coroniti fc Ka"t3 



1985 



2005 



1988'; Okuda' I2OO2I : lOkuda et al.ll2005l : lOhsuea et all 
Ohsuga 2007), and in SS 433 such large opening angles of 
the high-velocity flow were at odds with observations. The 
question in the previous studies was how to coUimate the 
high- velocity flow to a smaller angle. In the present study, 
this is not a problem. Indeed, the mass-outflow rates ob- 
tained for the relativistic axial outflow in the half opening 
angle of ~ 1° are sufficiently high to explain the observed 
mass-outflow rate of SS 433 jets, as far as the results of AD-1 
and AD-2 are considered. However, we have another ques- 
tion, why the relativistic axial outflow in our simulations is 
rather undistingishable from the broader high-velocity out- 
flow, while X-ray and optical jets of SS 433 are so brilliantly 
observed as the hot relativistic streams with a small open- 
ing angle in a region far beyond the present computational 
domain. 

We interpret the small opening angles of the jets in 
terms of the proposal (|Marshall. Canizares fc Schulzi[2003 : 
lFabrika|[200i ') that the observed expansion velocity in the 
transverse direction of the jets coincides with the sound 
velocity Ca in the region with a temperature ~ 10* K 
l|MarshaU. Canizares fc Schulj|2002l ) and that the half open- 
ing angle 9c observed in the distant region at r 3> 10^^ 
cm should be equal to ~ 2ca/Vjet, where Vjot = 0.26c. 
From our hydrodynamical results, in the axial flow we have 
Cs ~ 2.5 X 10~^c and T ~ 10* K, which agrees very well 
with the idea above. We suggest that, if future observations 
resolve the jets deeper to the central source, the jets (high- 
velocity outflows) will be viewed with large opening angles, 
> 1°. 

Studying H and He ab sorption l ines d uring the precess- 
ing motion of SS 433 disc, iFabrikal (|2004l ) derives the wind 
velocity Vw of SS 433 as a function of the angle from the 
disc plane: ~ (8000 ± 100) • sin^C + 150 km s"^ for ^ 
C ^ 30°. This gives velocities Vw = 400 - 2000 km s~^ (0.001 
- 0.007c) for 10° ^ C ^ 30°, which are smaller by a factor 
of 2 - 6 than the radial velocities K ~ 0.002 - 0.04c at the 
outer boundary in our simulations. The high wind velocity ~ 
1500 km s~^ (=0.005c) was sugges ted from analy sis of Hell 
emission near angles C ~ 70 - 80° (|Fabrikall20ol ). However, 
we obtain higher outflow velocities in these directions, 0.08 - 
0.1c. Therefore, we have to presume some mechanism of de- 
celerating the high-velocity flow in the distant region, such 
as interactions with the walls of gas cocoon surrounding the 
axial funnels where the jets expand. 

Although the absolute luminosity of SS 433 is not di- 
rectly observed, the kinematic luminosity L k is generally 
estimated to be 10 erg s" (|Kotanilll99^ ). which is con- 
sistent with ~ 1 - 2.5 xlO*" erg s~^ obtained in all cases. 
The absolute luminosity of SS 433 is interesting because 
it is the maximum luminosity among the accreting stellar- 
mass black holes. The compact objects in the recently ob- 
served ULXs are considered as stellar-mass black hole or 
intermediate mass black hole, depending on whether their 
luminosity far exceeds the Eddington luminosity Le or lies 
at sub-Eddington values (jMushotzkv, 2004 1 . The total lumi- 
nosities L obtained here are higher by a factor of 3 than 
that given by the estimate L/Le = 0.6 + 0.7 In mo ob- 
tained numerically in ID model l|Lipunoval Il999l ) and are 
higher by more than an order than the Eddington luminos- 



ity 1.5 X 10"^'' erg s~^. The coUimated outflows in our sim- 
ulations also lead to the outgoing radiation directed along 
the rotational axis. Actually, from the anisotropic radiation 
fleld obtained at the outer boundary for each model, we have 
the radial component Fr ~ 3 x 10^^ - 5 x 10^* erg s"^ cm"^ 
of the radiative flux F in the cone with a half opening an- 
gle 30° . If an observer within this cone assume the radiative 
flux to be isotropic over the entire surface, the apparent lu- 
minosity is 7 X 10"° - 10*^ erg s~\ Thus, SS 433 can be 
a representative of the supercritically accreting stellar-mass 
blac k hole cand i date, observed as the ULXs in nearby galax- 
ies llKing et all I2OOII: iFabrika. Abolmasov fc Karpovl I2OO6I : 



iBegelman. King fc Pringlell2006l ) 



We should pay attention to the remarkable modulations 
of A/edgc and in model AD-2. Their power spectra show 



QPOs-like signals at i/ ~ 4 x 10 



10"' and 0.5 - 2 Hz, 



but only the larger amplitude modulations survive as small 
modulations of the total luminosity with a quasi-period of ~ 
10 - 25 s, as mentioned in subsection 4.2. These instabilities 
develop as the recurrent hot blobs with variable size rising 
through the accretion disc and the polar funnel region. Such 
inhomogeneties could cause the luminosity modulations. Al- 
though occasional massive jet ejections, which are recog- 
nized as a clustering of flare events in radio light curves, 
were exhibited by Microq uasars, such as SS 433 (|Fiedleil 
Il987f ) and GRS 1915+105 (|Fosteij| 19961 ). X-ray observations 
of these events in SS 433 have been hardly performed so far, 
because the massive jet ejections are ra re, short, and ape- 
riodic. Nevertheless, iKotani et al.l (|2006l ) reported recently 
a variety of new phenomena, including a QPO-like feature 
near 0.1 Hz, rapid time variability, Oand shot-like activities, 
and suggested that an irregular massive jet ejection might 
be caused by the formation of small plasma bullets or knots 
in the continuously emanating flow. We propose that the ob- 
served QPO-like phenomena in SS 433 can be explained in 
terms of the recurrent hot blob phenomena found in model 
AD-2. 



7 CONCLUDING REMARKS 

We examined the jets and the disc of SS 433 with Mo ~ 
600Mc by time-dependent two-dimensional radiation hydro- 
dynamical calculations, assuming a-model for the viscosity. 
The initial discs are given by ID supercritical disc models 
with mass loss or advection. As the result, the total lumi- 
nosities obtained are 1 - 2.5 x 10*" erg s~^, which are 6 to 10 
times higher than the Eddington luminosity for M = WMq . 
From the initial models with advection, we obtain the total 
mass-outflow rates Mout ~ 4 x 10~^ and lO'^M© yr~^, and 
the relativistic axial outflow rates Mio ~ 1O~^M0 yr~^. 
These outflow rates agree well with the observed mass- 
outflow rates of the wind and the jets in SS 433. On the 
other hand, from the initial models with mass loss but with- 
out advection, we obtain the total mass-outflow rates and 
the ax;ial outflow rates smaller than or comparable to the ob- 
served rates of the wind and the jets respectively, depending 
on a. Still, while the mass-flow rate Mio of the axial outflow 
agrees well with the observed mass-outflow rate of SS 433 
jets covering the same half opening angle 1°, a problem is 
to be solved, why in the simulations the axial outflow is not 
distinguishable from the other high-velocity flow. 
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The calculated radial velocities of the high- velocity out- 
flow at the outer boundary around the disc are lar ger by a 
factor of 2 - 6 than the observed wind velocities ijFabrikal 
[2004 ). which are given as a function of the elevation angle 
from the disc plane of SS 433. Thus, for the radial veloci- 
ties to be consistent with the observed wind velocities, the 
outflow gas except the relativistic axial outflow must be de- 
celerated in the distant region beyond the present computa- 
tional domain by some mechanism, such as interaction with 
the inhomogeneous matter or the gas cocoon around the SS 
433 disc. 

The initial advective disc with large a evolves to the 
gas-pressure dominant, optically thin state in the inner re- 
gion. The inner disc generates instabilities in agreement with 
the previous stability analyses of the advection-dominated 
optically thin discs. As a result, we find remarkable mod- 
ulations of the mass-inflow rate at the inner edge and the 
disc luminosity. The modulations of Mcdgc present two types 
of variability; (1) the small amplitude variations with short 
time-scales of 0.5 -2 s, (2) the large amplitude variations 
with long time-scales of ~ 10 - 25 s. Only the variability (2) 
survives in the total luminosity because of the atmospheric 
absorption around the accretion disc. The disc instability 
results in the recurrent hot blobs, which develop outward 
and upward and produce QPOs-like phenomena in the total 
luminosity with the quasi-periods of ~ 10 - 25 s. The QPOs- 
like behavior of the luminosity and the hot blobs phenomena 
found here may explain the recent observations of a variety 
of new phenomena in SS 433, such as a QPO-like feature 
near 0.1 Hz, rapid time variability, and a shot-like activity 
ascribed to the formation of small plasma bullets. 
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